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Abstract 

We develop an empirical model capable of simulating the mean echo power cross prod- 
uct of CryoSat-2 SAR and SARIn mode waveforms over sea ice covered regions. The 
model simulations are used to show the importance of variations in the radar backscat- 
5 ter coefficient with incidence angle and surface roughness for the retrieval of surface 
elevation of both sea ice floes and leads. The numerical model is used to fit CryoSat-2 
waveforms to enable retrieval of surface elevation through the use of look-up tables and 
a bounded trust region Newton least squares fitting approach. The use of a model to fit 
returns from sea ice regions offers advantages over currently used threshold retracking 
10 methods which are here shown to be sensitive to the combined effect of bandwidth lim- 
ited range resolution and surface roughness variations. Laxon et al. (2013) have com- 
pared ice thickness results from CryoSat-2 and IceBridge, and found good agreement, 
however consistent assumptions about the snow depth and density of sea ice were 
not used in the comparisons. To address this issue, we directly compare ice freeboard 
is and thickness retrievals from the waveform fitting and threshold tracker methods of 
CryoSat-2 to Operation IceBridge data using a consistent set of parameterizations. For 
three IceBridge campaign periods from March 2011-2013, mean differences (CryoSat- 
2 - IceBridge) of 0.144 m and 1.351 m are respectively found between the freeboard 
and thickness retrievals using a 50% sea ice floe threshold retracker, while mean dif- 
20 ferences of 0.019 m and 0.182 m are found when using the waveform fitting method. 
This suggests the waveform fitting technique is capable of better reconciling the sea 
ice thickness data record from laser and radar altimetry data sets through the usage of 
consistent physical assumptions. 


1 Introduction 

25 Remote sensing records of Arctic sea ice thickness now span five decades and have 
shown nearly a two-fold decrease in mean winter thickness (Kwok and Rothrock, 2009), 
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while observations over the past three decades have shown a 1 7.2 % decade -1 decline 
in the areal coverage of multiyear ice (Comiso, 2012). The inter-related decline in sea 
ice thickness and multiyear ice coverage is tied to declining trends in ice age and sur- 
vivability (Maslanik et al., 2007; Maslanik et al., 2011). These changes have significant 
5 impacts to the climate, with a notable aspect of declining sea ice cover being linked to 
the observed higher than global average increase in Arctic surface air temperatures, 
a phenomenon known as Arctic amplification (Serreze et al., 2009). This occurs due 
to the increase of energy transferred from the atmosphere to the ocean as sea ice vol- 
ume decreases (Kurtz et al., 2011; Rigor et al., 2002), which enhances warming and 
io moistening of the lower troposphere (Boe et al., 2009; Screen et al., 2013). Changes 
in Arctic sea ice have also led to growing interest in determining predictability of the 
response of the sea ice cover to a changing climate. These interests range from ef- 
forts to improve short-term seasonal predictions (Lindsay et al., 2012; Eicken, 2013) 
to long-term predictions of when an ice-free summer may occur (Wang and Overland, 
is 2012), and if ice-free summers can be sustained over the long-term (Tietsche et al., 
2011). A key factor which links these disparate study areas is the need for continuous 
large-scale sea ice thickness observations to link physical processes to changes in sea 
ice and climate. 

The earliest historical remote sensing record of Arctic sea ice thickness is composed 
20 of declassified submarine sonar observations extending back to 1 958 (Rothrock et al., 
1999). The submarine sonar sea ice thickness record is composed of numerous pro- 
files within the central Arctic Ocean which need to be statistically analyzed to separate 
spatial, annual, and interannual variability within the limited regional data (Rothrock 
et al., 2008). Recent advances in satellite altimetry capabilities have enabled the de- 
25 duction of sea ice thickness and volume over the larger-scale Arctic Ocean basin on 
monthly time-scales extending from the beginning of the growth season in October to 
the beginning of the melt season in May. Laxon et al. (2013) produced the first results 
of Arctic sea ice thickness from ERS-1 and ERS-2 satellite radar altimetry measure- 
ments spanning October 1993 to March 2001 up to the latitudinal limit of 81.5°. The 
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ERS-1/2 radar altimetry record has also been extended using data from the Envisat 
satellite altimeter which showed large-scale thinning following the then record 2007 
sea ice minimum (Giles et al., 2008). For the period spanning 2003-2008, data from 
the ICESat satellite laser altimetry mission provided a record of sea ice volume with 
5 increased coverage up to the latitudinal limit of 86°. ICESat observed further decline 
in the thickness and volume of the Arctic sea ice cover in agreement with the radar 
altimetry results (Kwok et al., 2009). Presently, ESA’s CryoSat-2 mission (Wingham 
et al., 2006), launched in 201 0, is producing a continuous time series of radar altimeter 
measurements up to a latitudinal limit of 88°, providing unparalleled coverage of the 
io Arctic sea ice cover. 

Laxon et al. (2013) produced the first estimates of sea ice thickness and volume 
derived from CryoSat-2 data and validated the data with multiple in-situ data sets. The 
CryoSat-2 results were combined with ICESat estimates to produce the first decadal- 
scale record of basin-wide Arctic sea ice volume from satellite altimetry. Data were 
is also compared to estimates from the Pan-Arctic Ice-Ocean Modeling and Assimilation 
System (PIOMAS) model, which has shown volume loss of nearly 3x 10 3 km 3 decade -1 
from 1 979-201 0 (Schweiger et al., 201 1 ) and similar trends from 1 979 to the present. 
The combined ICESat and CryoSat-2 time series of sea ice volume change provides 
a useful tool to assess the PIOMAS data set which shows a loss of sea ice volume 
20 over a much longer time period. 

With the advent of sea ice volume records from different satellite altimetry data 
sources comes the need to reconcile the assumptions used in the retrieval processes 
to produce a continuous time series and quantify uncertainties. Differences in sea ice 
thickness estimates from altimetry data arise in particular to the use of different density 
25 values and snow depth estimates which are used in the retrieval of sea ice thickness. 
These quantities are due to environmental processes and should be applied in a con- 
sistent manner in the retrieval of sea ice thickness regardless of which instrument is 
used. In the case of sea ice density, previous studies have utilized a wide range of 
values, which will result in large differences between data sets if the same physical 
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assumptions are used. For example, in the study by Kwok et al. (2009) an ice-density 
of 925 kg m -3 was used, Kurtz et al. (2011) used a value of 91 5 kg m -3 , while Laxon 
et al. (2013) used an estimate of 91 7 kg m -3 for first year ice and 882 kg m -3 for mul- 
tiyear ice. In these studies, the range of sea ice density values for multiyear ice is 
5 particularly large at 43kgrrT 3 . For a typical multiyear sea ice floe with 60 cm of snow- 
ice freeboard and 35 cm of snow depth, the sea ice thickness estimate differs by 1.1 m 
within this range of ice densities. Despite the large-scale mean agreement of the sea 
ice thickness data sets described in previous studies, this discrepancy in physical as- 
sumptions points to the source of the differences as being due to potential biases in 
io the freeboard and snow depth data sets used. This large discrepancy underscores the 
need to establish a set of consistent physical constants for use in the retrieval of sea 
ice thickness from satellite radar and laser altimetry data. 

The focus of this study is to develop a new method for the retrieval of sea ice free- 
board from CryoSat-2 data. We demonstrate that this method is consistent with inde- 
15 pendent measurements from airborne laser and radar altimetry data sets from NASA’s 
Operation IceBridge mission to retrieve sea ice thickness which eliminates the need to 
utilize different ice density and snow depth values as an effective bias correction. The 
study is organized as follows: Sect. 2 describes the data sets used in the study. Sec- 
tion 3 describes the physical model which is used to simulate CryoSat-2 returns from 
20 Arctic sea ice covered regions. A procedure to utilize the model to fit CryoSat-2 wave- 
forms for the retrieval of surface elevation is developed in Sect. 4. The new retrieval 
procedure is used to estimate Arctic sea ice freeboard and the results are compared 
to a threshold tracking method and independent freeboard observations from airborne 
data in Sect. 5. The results are then summarized and future improvements to the re- 
25 trieval method are described in Sect. 6. 
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2 Data sets 

The primary data set used in this study comes from ESA’s CryoSat-2 satellite (Wing- 
ham et al., 2006). Data are taken from the baseline B Level IB SAR and SARIn 
mode data products for March 2011, 2012, and 2013. Example CryoSat-2 SAR mode 
5 waveforms and terminology employed in the description of the waveform features are 
shown in Fig. 1. CryoSat-2 is a radar altimeter which operates at a center frequency 
of 13.575 GHz and has a receive bandwidth of 320 MHz. The effective footprint size af- 
ter post-processing is ~ 1650 m in the across track direction and ~ 380 m in the along 
track direction. The power detected echoes contain 128 range bins in SAR mode and 
io 512 range bins in SARIn mode, each range bin is sampled at 1.563 ns (0.234 m range 
resolution in vacuo). The satellite operates in SARIn mode over a spatially limited sec- 
tion of the Arctic Ocean, SARIn mode utilizes dual receive antennas to obtain phase 
information which can be used to detect the angle of off-nadir reflections. The focus of 
this study is to describe retrieval methods which can be used for the power detected 
is waveforms. Thus, in order to maintain consistency in the retrieval algorithms developed 
here, phase information is not used and the SARIn data are truncated from 512 to 128 
range bins. 

The window delay field in the level IB data provides the one-way travel time from 
the center range gate to the satellite’s center of mass. We use this data to retrieve an 
20 elevation above the WGS84 ellipsoid by multiplying by the speed of light in vacuum, 
applying geophysical and retracking corrections, and subtracting these from the satel- 
lite center of mass altitude. Corrections for the elevation are given in the data products 
for the wet and dry tropospheric delay time, ionospheric delay, oscillator drift, inverse 
barometer effect, dynamic atmospheric correction, ocean equilibrium tide, long period 
25 ocean tide, load tide, solid earth tide, and pole tide. These corrections have been ap- 
plied to the data used in this study. Retracking the mean scattering surface within the 
radar waveform is the focus of this study which allows for the sea ice surface elevation 
to be determined. 
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Data from NASA’s Operation IceBridge airborne mission are used for com- 
parison with monthly mean CryoSat-2 data for three campaign periods spanning 
March 2011 to March 2013. We use data from the IceBridge Sea Ice Freeboard, 
Snow Depth, and Thickness products (Kurtz et al., 2012) from 16-28 March 2011, 
and 14 March-10 April 2012. Data from the quick look data set have been used for 
20-27 March 2013 since final data from the campaign are not yet available. The 2013 
quick look data utilizes new processing techniques to minimize freeboard biases (Kurtz, 
2013), but it is possible that the uncertainties in this data set are higher than in the 
final archival product (Kurtz et al., 2013a). The data consist of measurements from 
the Airborne Topographic Mapper laser altimeter (Krabill, 2010), Digital Mapping Sys- 
tem camera (Dominguez, 2010), and the University of Kansas’ 2-8 GHz snow radar 
(Leuschen, 201 0; Panzer et al., 201 3). Data from the individual instruments have been 
synthesized to provide sea ice freeboard, thickness, and snow depth at a 40 m spatial 
sampling resolution along all available flight lines using the methodology described in 
Kurtz et al. (2013b). Uncertainty estimates are also provided with the data products, 
which are estimated from the number of sea surface height tie points, distance to the 
local sea surface tie points, and the estimated covariance of the sea surface height for 
each flight. In this study we restrict data usage to where the uncertainty in the laser 
altimeter derived sea ice freeboard is less than 0.1 m. A map of the IceBridge derived 
ice freeboard used in the study is shown in Fig. 2. 


3 CryoSat-2 multi-look echo phenomenology 

In this section the behavior of the CryoSat-2 waveforms over surface types encountered 
in sea ice covered regions of the Arctic are simulated through the use of a physical 
model. The model shows the theoretical variation of the echo tracking point needed for 
the retrieval of surface elevation from the different surface types encountered. Before 
describing the model used in the simulation of CryoSat-2 returns, we acknowledge that 
due to the inherent complexity of scattering from sea ice covered regions assumptions 
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need to be made to simplify the problem to attain a tractable solution. In particular, 
we treat the scattering from sea ice as a surface problem. We furthermore assume the 
height deviations within the radar footprint are Gaussian and have a spatial exponential 
autocorrelation. Where appropriate, we note in the text where assumptions have been 
5 made and attempt to justify them, though there is of course a limit to the accuracy of 
the assumptions used. Given the assumptions made in the scattering model, it must 
still be treated as empirical, the validity of which is thus based on the degree to which 
it is able to model the phenomena of CryoSat-2 returns. Towards this end, potential 
improvements to the physics of the model are discussed in Sect. 6. 

10 3.1 Physical model for CryoSat-2 echoes 

Here we provide the theoretical basis for modeling the mean echo power cross product 
from CryoSat-2 SAR and SARIn mode waveforms over sea ice covered regions of the 
Arctic. CryoSat-2 differs from previous generation pulse-limited radar altimeters (e.g. 
Envisat, Jason-1/2) largely in two ways: (1) the radar altimeter of CryoSat-2 consists of 
is two antennas which have been narrowed in the across-track dimension to fit within the 
launcher fairing, thus it has an elliptical rather than circular antenna pattern which alters 
the impulse response (Wingham and Wallis, 2010). (2) Unfocussed aperture synthesis 
is employed to reduce the along-track footprint size of the surface return. The level 1b 
data products are the result of a beam formation process which sums phase weighted 
20 and slant range corrected echoes taken from different look angles (see Wingham et al., 
2006 for details). 

The received radar echo, H'(t), from a uniformly backscattering planar surface can 
be expressed as (e.g. Brown, 1977; Raney, 1998; Wingham et al., 2004) the convolu- 
tion of the compressed transmit pulse after signal processing, P t (r), the surface height 
25 probability density function, p(r), and the “rough surface” impulse response (a factor of 
the surface geometry and antenna pattern), /(r), 
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where ® represents convolution and r represents the time delay relative to the time of 
the first surface arrival (r = 0 is thus the point of closest arrival, which is here consid- 
ered to correspond to the mean scattering surface). The use of Eq. (1) assumes that 
only surface scattering from the snow-ice interface is present (i.e. the surface is as- 
sumed to be perfectly conducting), surface scattering from the snow-air interface and 
volume scattering from within the snow and ice layers are neglected. This assumption 
is justified when the dominant reflection of energy occurs from the snow-ice interface 
i.e. when the density of the snow pack is relatively low (compared to ice), the reflection 
coefficient of the snow-air interface is much less than that of the snow-ice interface, 
and volume scattering within the ice layer is low. This will occur in practice when the 
snow pack does not contain wet or icy layers, and the roughnesses are not sufficiently 
different so as to cause a significant increase in the ratio of the Fresnel reflectivities of 
the snow-air and snow-ice interfaces. 

The CryoSat-2 compressed transmit pulse is well-represented by a sine function 
described in Galin et al. (2013) as 

PM = P 0 sinc 2 (7r£ w r) (2) 

where p 0 is the peak power of the compressed pulse and B w is the received bandwidth 
(£ w = 320 MHz). 

The surface height probability density function is not known a priori, rather it must be 
determined through analysis of the waveform shape as will be shown in Sect. 4. Here 
we assume that p(r) follows a Gaussian distribution given by 

p<T)= ^ exp (‘H£0 <3) 

where o c = 2o /c which is the surface roughness in the time domain and c is the speed 
of light in vacuo. A Gaussian height distribution was chosen since it is dependent on 
only a single parameter, and it is not presently known what form the height distribution 
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of the snow-ice interface will take over the radar footprint. Future research in this area 
is needed to determine whether a different height distribution assumption can be used 
to improve the accuracy of the retrievals. 

Following Brown, (1977), the impulse response used to determine the power for 
a conventional altimeter is 


Ut) = 


64n 3 


/, 


Illuminated area 


g{0,Cl)o°(r,x)5 




( 4 ) 


where X is the center wavelength, g is the antenna pattern, cr° is the radar backscat- 
tering coefficient, t is the time from the instant of transmission, r is the range from the 
radar to the elemental scattering area 6A on the surface. The angular component of 
the antenna gain pattern, (0,fi), is measured relative to the antenna boresight, while 
the angular component of the backscatter dependence, (T,/) is relative to the surface 
normal. To simplify the problem, we here make the approximation that the satellite pitch 
and roll are zero and that the surface normal is parallel to the nadir direction. This al- 
lows for the angular components to be written as (0,fi) = (T,/) = (y,®) where y and 
® are respectively the polar and azimuth angles subtended at the altimeter between 
the antenna boresight and scattering element. The standard deviations of the pitch and 
roll values over Arctic sea ice regions are small at 0.006° and 0.01°, respectively, with 
some of the observed variability due to noise in the star tracker measurements (Galin 
et al., 2013). The recorded mean pitch and roll over Arctic sea ice regions is less than 
0.01°, however, there is a known bias in the recorded pitch and roll values due to an 
error in the star tracker rotation matrices (Galin et al., 2013; Galin et al., 2014) which 
should be taken into account if a more physically exact characterization of the impulse 
response is desired. 

Following Wingham et al. (2006), Eq. (4) can be extended for application with 
CryoSat-2 through the addition of a synthetic beam gain term, where d 0 {cos£ -sin^) 
is the synthetic beam gain, which is a function of the angle between the direction of 
a scattering element and the satellite velocity vector, £, and the look angle of synthetic 
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beam k from nadir, £ k . As described in Galin et al. (2013), the impulse response must 
also be summed over the different look angles used in the beam formation process. 
For the study by Galin et al. (2013) the look angle was defined in terms of the higher 
angular sampling in the burst which allows for the retrieval of the echo power and thus 
5 a 0 . Since we are only concerned with the echo power shape, in the context of this 
study we define £ k = k- 0.0238 which refers to the look angles from the stack data (the 
definitions of “burst” and “stack” are described in Wingham et al., 2006). Similarly, d 0 is 
the FFT of a Flamming window which is used in the formation of the mean echo cross 
product and D 0 is the FFT gain of the synthetic aperture minus the Flamming window 
io loss. Using these definitions, the impulse response can then be written as 


l s (t) = 


Do£ 

647T 3 


(/V b -1)/2 

I 

/f=-(/V b -1)/2 


Illuminated area 


g{Y,GD)o°{Y,GD)S 




(cos f — sin £ k )-—dA (5) 


We expand on the models of the CryoSat-2 impulse response described in Wingham 
et al. (2006) and Galin et al. (2013) by including a backscatter coefficient which varies 
is with incidence angle, this will be shown to be necessary for modeling of the CryoSat-2 
waveforms over sea ice. Over the range of incidence angles (up to 0.76°) encountered 
by CryoSat-2 for sea ice regions, the type of scattering is here assumed to remain 
within the specular scattering regime. For specular scattering, only surface facets which 
are tilted normal to the direction of the incident radiation contribute to the backscatter- 
20 ing (Flagfors, 1964; Valenzuela, 1977). Hagfors, (1964) showed that for smoothly un- 
dulating surfaces with an exponential autocorrelation of height features the theoretical 
received power, O, with respect to incidence angle, 0 , for a plane wave undergoing 
specular scattering is 


CD (0) = cp(0°) 


COS 4 0 + 


,2 k 0 h 2 m 


2 

sin 2 0 


^3 

2 


( 6 ) 
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where h m is the rms height deviation, / is its length scale, and k Q is the carrier 
wavenumber. We relate this to variations in a 0 (y,®) by considering a scenario which 
varies only in the incidence angle by taking a ratio of the received power for nadir 
scattering and at an angle <p = y (which considers only variations in the polar angle 
direction) which gives 


cb(0°) _ <j(0°) 
W) ” o{Y) 


cb(0°) 


<t>(0°) 


cos 4 y + 


(2/^J 


sin 2 y 


( 7 ) 


For the small incidence angles encountered by CryoSat-2, we assume cos 4 y ~ 1 which 
gives the approximate variation of backscatter with incidence angle to be 

-3 

o°(y) = o°(0°) [l +a'sin 2 yj 2 (8) 

where a°(0°) is the backscattering coefficient at nadir and a = ( — is a dimension- 

less variable that quantifies the efficiency of backscattering from a surface as a func- 
tion of incidence angle, a is not known a priori and is determined from estimates of the 
waveform shape as described in Sect. 4. Note this treatment does not quantify o°(0°) 
and its dependence on the rms height deviation and length scale, it is here treated as 
an unknown constant which affects the echo amplitude, but not shape. Thus, the value 
of o°(0°) is not retrieved in the context of this study. In choosing Eq. (8) a Gaussian 
height distribution with an exponential autocorrelation of height features is assumed 
since it allows for an explicit relationship between the surface height deviation and 
surface slope to be realized (Hagfors, 1964), which then allows for the simple mathe- 
matical characterization of the power directionality dependence described in Eq. (6). 
Physically, we may expect this assumption to be valid when the surface consists of 
a largely homogenous and isotropic field of height features. However, in areas such 
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as heavily ridged ice this assumption will likely introduce additional uncertainty in the 
results since the height distribution may not be Gaussian and the autocorrelation may 
have a different spatial dependence. 

Equation (5) can be reduced to a line integral around an isorange circle following 
the approximations to the scattering and geometry described in Wingham et al. (2006) 
and Wingham, (1 995). The expression used here for the impulse response of the multi- 
looked echo power shape follows the expression described in Galin et al. (2013) with 
the additional use of Eq. (8) to include the variation of backscatter with incidence angle. 
Consequently, the full expression for the impulse response is written as 


X 2 GqD 0 co°(0°) (/Vb _l )/2 / qhfr 

/ ( T )= — J H It + —— 


32n 2 h 3 q 


( 9 ) 


k=-(N b - 1)/2 


exp 


-2 A 


_/± + ±VjW T+ ^ N1 

n 2 Vk? y 2 J h Y 2 \ c , 


2n 

I 


d0exp 




\| hq 3 


T + 


nh£l\ /l 1 \ 2cos(2 9) f nhff 


cos 9 ( — + — 1 - 


y i yI 


hqY 2 2 


T + 


/ 


1 + 


V 


a 

h 2 


/ '^) 2 + |^ t+ ^ +2 


H k 


a ch l nh£ 

N n V c 


2x\^- 3/2 


\ 


/ 


2 ^0.54 - 0.46 cos ~ ) cos ( 2 Vs("- ^k COs6 ~ 


where Table 1 summarizes the parameters and symbols used in the equation. The ex- 
pression after the summation in the second and third lines of Eq. (9) correspond to 
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the contribution of the elliptical antenna pattern, it is taken from Wingham and Wal- 
lis (2010). The fourth line of Eq. (9) corresponds to the variation of backscatter with 
incidence angle. The fifth line corresponds to the gain of each synthetic beam and 
the application of a Hamming window. The form of the equation accounts for the slant 
range correction of each synthetic beam which is employed in the CryoSat level IB 
data processor. 

3.2 Waveform simulations 

Equations (1), (2), (3), and (9) describe the physical model used to simulate CryoSat-2 
waveforms over sea ice. As there is no closed form solution to these equations they 
must be calculated numerically. Here we detail the theoretical behavior of the waveform 
shape over both sea ice leads and floes by using this model. The free parameters in 
the model simulations shown here are o and a. o is varied from 0 to 0.4 m which rep- 
resents the expected range from a smooth lead to ridged sea ice over the CryoSat-2 
footprint, a is varied from 0 to 5 x 10 7 which represents the range from open ocean re- 
turns and very rough sea ice (where a « 0; there is little to no backscatter dependence 
with incidence angle) to a perfectly smooth lead where backscatter from the nadir point 
dominates the echo. The leading edge of the waveform is affected by both a and the 
surface roughness. Increasing surface roughness increases the width of the echo, par- 
ticularly from the rise time edge to the peak, it has much less impact on the trailing 
edge (Wingham et al., 2004) which is largely affected by the r 2 behavior of the area of 
the range cells (Wingham et al., 2006). Equation (9) shows that as a becomes large, it 
dominates the decay of the trailing edge of the waveform. 

3.2.1 Sea ice leads 

The effect of variation in o° with incidence angle is to decrease the effective illuminated 
area on the surface. This is apparent over sea ice leads where returns from geometri- 
cally small leads dominate the echo from radar altimeters (Drinkwater, 1991). Figure 3 
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shows simulations of CryoSat-2 waveforms for o = 0.02 m and the observed range of 
a over sea ice leads (shown in Sect. 5). It can be seen that placement of the tracking 
point to determine the surface elevation for sea ice leads is sensitive to a, because 
it determines the contribution of the off-nadir beams used in the retrieval. The mean 
scattering surface can be seen to correspond to the echo maxima for a > 5 x 1 0 7 , and 
progressively moves toward the waveform leading edge as a decreases. 

Mathematically, it can be shown that over smooth leads {o ~ 0; p(r) = <5(r)) with 
a suitably large value of a the impulse response function goes to a delta function, /(r) = 
<5(t, £ k ). Note however, that Eq. (6) assumes incoherent reflections whereas scattering 
from a surface with o = 0 will be coherent, which will affect the pulse amplitude, but not 
shape, which is the focus of this study. The received echo shape for a perfectly smooth 
lead will be 

Y(T) = P f (T)®,5(T,^)®<5(T) = P,(T) 

which is simply a copy of the transmit pulse. This is also illustrated in the lowest am- 
plitude waveform in Fig. 3 and can be seen in select CryoSat-2 waveforms over leads 
(an example of which is shown in Sect. 4). As a decreases, returns from off-nadir are 
incorporated and the trailing edge of the waveform becomes longer. With the inclu- 
sion of more returns from off-nadir the mean scattering surface shifts leftward from 
the maximum peak power. Quantitatively, for a = 0.02 m the echo peak corresponds to 
r = 0.000 ns for a = 5 x 1 0 7 , and r = 0.203 ns (0.030 m) for a = 5 x 1 0 5 , this range of 
3 cm is the maximum sensitivity of elevation retrievals from leads due to a variations. 

3.2.2 Sea ice floes 

Simulated CryoSat-2 echoes from sea ice floes are shown in Fig. 4. For a Gaussian 
surface height distribution, Fig. 4 shows the mean scattering surface occurs when the 
leading edge reaches approximately ~ 85-95 % of the peak value, with some variation 
of this threshold due to a and o variations. The result that the retracking point for SAR 
echoes is near the peak, rather than at the half power point as is found in conventional 
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pulse-limited altimeters, was shown previously by Wingham et al. (2004). One point to 
note is the effective 1 .563 ns sampling resolution of the instrument may not allow for the 
peak power to be well-determined for waveforms with low o and high a, this will impact 
threshold algorithms which rely on a peak power ratio. For the 50% threshold tracker 
used by Laxon et al. (2013), the simulations show biases due largely to variations in o 
and less significant biases due to a variations. The biases range from -2.969 ns (lead- 
ing to an elevation bias of +0.445 m) for a - 0.4 m and a = 10 3 to -0.531 ns (+0.08 m) 
for a = 0 and a = 10 5 . The variation of the r = 0 point for different threshold values 
shown in Fig. 4b demonstrates that the freeboard for threshold tracking methods will 
likely be biased. However, the basin-wide bias encountered in an operational setting 
can not be accurately quantified from the simulations since it will be dependent on 
the combined effect of the surface roughness, the surface height distribution within the 
footprint, and the finite range resolution of the instrument. 


4 Surface elevation retrieval algorithm 

In this section, the physical model is combined with a least squares fitting procedure 
to estimate the mean scattering point and mean surface roughness within CryoSat-2 
echoes from varying surface types. This least squares fitting procedure is analogous to 
routines which fit physical models to waveforms over ocean returns to retrieve surface 
elevation and other parameters such as significant wave height. Since Eq. (1) as de- 
veloped here does not have a closed form solution, we describe the procedures which 
are used to fit the waveforms from numerical solutions. We show that through the use 
of look-up tables, the computation time of a least squares fitting routine is sufficient to 
fit the waveforms without the need for a closed form solution. Our fitting routine can fit 
a single CryoSat-2 LI B SAR/SARIn waveform on the order of one to ten seconds using 
a standard desktop computer, and a single month of CryoSat-2 data over Arctic sea 
ice can be processed in ~ 10 days. Thus, the retracking method using the best model 
fit is practical from a processing standpoint. 
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4.1 Fitting routine 

In order to speed up calculation and enable fitting of individual waveforms, we calculate 
a look-up table of L{r) = P t (r) ® /(r) for a discrete set of cases encountered over Arctic 
sea ice and placed the data on an irregular grid. A fitting routine with pre-computed 
interpolation coefficients is then used to linearly interpolate between these discrete 
cases and quickly provide a solution for the function and its first and second deriva- 
tives (using the method of finite differences) for any queried point within the parameter 
space. 

After creation of the look-up table, a least squares fitting routine using a bounded 
trust region Newton method (MATLAB function Isqcurvefit; described in Coleman and 
Li, 1996) is used to minimize the difference between the model fit and each CryoSat-2 
echo power waveform, P r A bounded trust region Newton method was chosen because 
the method is globally convergent, relatively independent of the problem size, and few 
iterations are needed to converge to a solution (Coleman and Li, 1994; Coleman and 
Li, 1996). Four free parameters are used in the fitting routine which is characterized by 
the equations 

P m (r) = A f L(T,a,o)®p(T,o) (10) 

128 

min J [ P m (T/) - P r (T, + t)f (1 1 ) 

/'= 1 

where P m is the modeled waveform, P r is the observed echo power, and r, corresponds 
to the observed echo power at point / of the waveform. The four free parameters are: 
(1) A f , the amplitude scale factor, (2) t the echo time shift factor, (3) a, (4) o. 

Given the dynamic range of the input parameters, and the fact that the solution which 
minimizes the square of the differences may not be physically correct, we specify an 
initial guess for each waveform and provide upper and lower bounds for the unknown 
parameters which are dictated by the physical system. The initial guess for A f is taken 
to be equal to the waveform peak power for all cases. For all other parameters, the 
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methods for initial guess and upper and lower bounds are provided in the specific 
cases outlined below. 

C/) 

GO 

4.2 Leads 

TJ 

0) 

CryoSat-2 data over leads are identified in a similar manner to Laxon et al. (2013). 

First, a pulse peakiness parameter is calculated as 
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128 

PP= Jma x(P t )/P t (i) (12) 

/= i 

Following Laxon et al. (2013), leads are defined as having a PP > 0.18 and a stack 
standard deviation < 4. An initial guess of o = 0.02 m is used and the bounds are taken 
to be 0 < a < 0.1 m. f, is first estimated to be the point of maximum power. An initial 
guess for a is estimated from the theoretical waveform peak to tail ratio which is taken 
from Sect. 3 and shown in Fig. 5. The tail is defined as the mean power of the six 
range bins (1 0 ns) following the point of peak power. The bounds for a are taken to be 
^ < a 0 < 100a 0 where a 0 is the initial guess. 

Example fits to CryoSat-2 waveforms over leads are shown in Fig. 6. The figure 
shows the behavior of the CryoSat-2 waveform for increasing values of a. As shown in 
Sect. 3.2.1 , over smooth leads (a ~ 0) with a large value for a (a £ 5x10 7 ), the received 
waveform is simply a copy of the transmit pulse which may be slightly broadened by 
the small surface roughness within the lead. An example CryoSat-2 waveform showing 
this behavior can be seen in Fig. 6d. For all lead cases, the tracking point for the mean 
scattering surface is near the maximum peak of the return. For a < 5 x 1 0 7 , returns from 
off-nadir begin to broaden the waveform and shift the mean scattering surface leftward 
from the maximum peak due to the inclusion of off-nadir look angles i.e. it is determined 
by the combined effect of the impulse response for each off-nadir look angle and the 
slant range correction used in the data processor. As shown in Sect. 3.1.1, the tracking 
point for the mean scattering surface over leads is thus sensitive to the choice of a 
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with a maximum uncertainty of 3 cm, but the uncertainty in surface elevation caused 
by errors in the choice of a in the fitting routing are likely small since the returns from 
leads can be seen to be very well-represented by the physical model. 

4.3 Sea ice floes 

For sea ice floes we follow Laxon et al. (2013) and define waveforms from floes as 
having a PP < 0.09 and a stack standard deviation greater than 4 (3 for SARIn mode 
regions). The initial guess for t, is taken from Laxon et al. (2013) as the first point 
where the waveform power reaches 50% of the power of the first peak. The first peak 
is defined as the first local maximum on the waveform leading edge with a power value 
greater than 50% of the point of highest power in the waveform. The waveform is only 
used when the power of the first peak is greater than 80% of the highest power value 
in the waveform. The upper and lower bounds for t, are taken to be ±6 ns (±0.9 m) 
from the initial guess point. The initial guess for a is determined by the trailing edge of 
the waveform in a similar manner to that of leads. Figure 7 shows the tail to peak ratio 
which is used for sea ice floes. For sea ice floes, the tail is taken from the mean power 
of the set of measurements between 90 to 120 ns (58-78 range bins) after the point of 
peak power. The upper bound for a is taken to be 1 00 times the initial guess for a and 
the lower bound is taken to be the initial guess for a divided by 100. The trailing edge 
of the waveform is used for sea ice floes since the larger off-nadir angles experienced 
at larger delay times tends to eliminate the more “peaky” aspects from flat targets such 
as off-nadir leads. The initial guess for o is set to 0.1 m, with a range of possible values 
from 0 to 1 m, if the initial guess for a is less than 8000 (which occurs over the open 
ocean and very rough sea ice floes) then the upper bound for o is set to 6 m. 

Example fits of the physical model to sea ice floes are shown in Fig. 8. The model fits 
the CryoSat-2 return waveform very well for both smooth (Fig. 8a) where o = 0.05 m) 
and rough ice (Fig. 8b where a = 0.34 m), which provides confidence in the ability of fit- 
ting model to be used to retrieve surface elevation over sea ice floes. We note that while 
Sect. 3 estimated the tracking point to be where the waveform leading edge reached 
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~ 80 % of the peak power, the fitted results shown in Fig. 8 demonstrate the finite range 
resolution of the instrument changes this value in practice such that a choice of an 80 % 
threshold would not yield accurate results in all cases. This is illustrated in Fig. 4b for 
the a - 0 cases where the waveform leading edge, peak power point, and trailing edge 
are all located within the effective 1.563 ns sampling resolution of CryoSat-2, thus the 
peak power may not be adequately determined within a given waveform due to sam- 
pling limitations. 

Model fits for areas with a mixture of smooth and rough surface types are shown 
in Figs. 8c and d. The physical model developed in Sect. 3 assumed a surface with 
uniform characteristics which leads to the observed variations from the model fit. When 
the surface is not largely homogeneous within the CryoSat-2 footprint, a mixed return 
will result due to the different backscattering properties within the footprint. This is due 
to the inter-related variations in o°, a, and o which will combine to create a signal 
which has multiple peaks, unlike the single-peak smooth theoretical echoes seen in 
the model. The use of the pulse peakiness parameter to distinguish between sea ice 
floes and leads is discussed in Peacock and Laxon (2004), and Laxon (1994). In this 
study, we used the pulse peakiness and stack standard deviation thresholds used by 
Laxon et al. (2013) to minimize errors caused by mixed returns. The fitted returns show 
that if the smooth areas within the radar footprint have a large enough off-nadir angle 
so as to make the secondary peaks distinguishable from the main peak, then they do 
not largely impact the fitting routine since the location of the mean scattering surface is 
on the waveform leading edge. 


5 CryoSat-2 derived sea ice properties 

In this section we discuss the procedure for retrieving sea ice properties, including 
freeboard, roughness, and thickness from the CryoSat-2 data set. The freeboard and 
thickness results are compared to a threshold tracker method for sea ice floes and to 
independent measurements from NASA’s Operation IceBridge campaign. 
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5.1 Sea ice property retrievals 

O 

C 

For the frequency range used by CryoSat-2, the surface return from sea ice covered 
regions is assumed to be from the snow-ice interface as has been shown to be the 
dominant reflecting surface in laboratory experiments (Beaven et al., 1995), and field 
experiments when cold, dry snow is present (Willatt et al., 2011). Sea ice freeboard is 
thus defined here to be the height of the ice layer a.s.l. and is calculated as 
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CryoSat-2 sea ice 
freeboard 


f^radar — ^fioe ^ssh 0 3) 

where h ftoe is the sea ice floe elevation and /7 ssh is the sea surface elevation. For all 
CryoSat-2 waveforms, we first removed the time varying sea surface height parame- 
ters outlined in Sect. 2. We then apply the retracking correction which is taken from 
the waveform fitting model used in Sect. 4. We calculate a monthly mean freeboard 
by gridding all sea ice floe and lead data points to a 25 km polar stereographic grid, 
with each grid point required to contain five floe elevations and five sea ice lead eleva- 
tions to be flagged as containing a valid freeboard retrieval. After this initial gridding we 
then smooth the data by taking the average value for all points within ±2 grid points. 
This effectively reduces the spatial resolution to ~ 1 00 km. A map of the mean gridded 
CryoSat-2 freeboard retrievals is shown in Fig. 9. Since the radar measures the ice 
freeboard which is the dominant factor in the retrieval of sea ice thickness, the spatial 
distribution of freeboard heights is expected to be similar to that of the ice thickness. 
The map shows a spatial pattern which is consistent with past observations (Bourke 
and Garret, 1987; Kurtz et al., 2011; Kwok et al., 2009; Laxon et al., 2013) with high 
freeboards in the multiyear ice regions north of the Canadian Archipelago and Green- 
land, and lower freeboards in the first year ice regions of the Arctic Ocean and outlying 
seas. 

The roughness of the scattering surface, o, can also be retrieved from the waveform 
fitting method. A map of the surface roughness (excluding sea ice leads) is shown in 
Fig. 10. The sea ice floe roughness also corresponds well to what may be expected 
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from known dynamics and circulation patterns in the Arctic Ocean with the rough- 
est ice corresponding to the multiyear ice area north of Greenland and the Canadian 
Archipelago. Gridded data points for o and logger are highly correlated with a correla- 
tion coefficient of -0.8, this demonstrates that, as expected, an increasing ice surface 
roughness corresponds to a lower angular variation in the radar backscatter coefficient. 

In order to retrieve the ice freeboard needed for sea ice thickness retrievals, a geo- 
physical correction to the CryoSat-2 freeboard must also added to account for vari- 
ations of the speed of light within the snow pack on sea ice. This is given as fb = 
fbradar + h c , where the correction factor, h c , is given as 


h c = h s 



(14) 


where h s is the snow depth and c snow is the speed of light within the snow pack. c snow 
is parameterized following Tiuri et al. (1984) to be 


\/l + 1.7p s + 0.7p 


(15) 


where p s is the density of snow with units of gem -3 . For the comparison with IceBridge 
data discussed in the next section, this geophysical correction adds a mean value of 
4.9 cm to fb radar to attain the true ice freeboard. 

Sea ice thickness, h v can be retrieved from the CryoSat-2 data set through the as- 
sumption of hydrostatic balance 


= 


Pw 

Pw — Pi 


fb + 



(16) 


where p w , p v p s are the respective densities of sea water, ice, and snow. Thus, the 
retrieval of sea ice thickness requires an independent snow depth data set as well 
as assumptions of the density properties of the surface. In this study we use density 
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assumptions which are discussed in Kurtz et al. (2013b) to be consistent with the Ice- 
Bridge data. The density of sea ice is taken to be 91 5 kg rrf 3 , the density of snow is 
taken to be 320 kg m -3 , and the density of sea water is taken to be 1024 kg m -3 . Using 
these values Eq. (16) can be written as 

s h t = 9.39fb + 2.94h s (17) 

5.2 Comparison of CryoSat-2 freeboard and thickness data 

In this section we compare CryoSat-2 retrieved freeboard data using the waveform 
fitting method and an empirical lead and threshold floe tracker. We then provide an 
io independent comparison to Operation IceBridge data. 

5.2.1 Comparison to sea ice floe threshold tracker 

In order to illustrate differences between the new freeboard retrieval method developed 
in this study, we compare freeboard retrievals from the waveform fitting retracker to 
a similar freeboard retrieval method outlined in Laxon et al. (2013). The method of 
is Laxon et al. (2013) uses a 50% threshold tracker in the retrieval of sea ice floe el- 
evations, and retracks sea ice lead returns using an empirical fit function described 
in Giles et al. (2007). Our reproduction of a similar method is hereinafter referred to 
as the ELTF (empirical lead and threshold floe) retracker. We note that several differ- 
ences are present between the freeboard retrieval used by Laxon et al. (2013) and the 
20 ELTF method used here. The primary difference is that Laxon et al. (2013) subtracted 
a bias from the sea ice lead elevations by taking the difference between returns from 
the ocean when sea ice is not present and returns from leads in the nearby ice pack. 
This was done following Giles et al. (2007), but was not done in the ELTF freeboard 
retrievals. Additional differences include (but are not limited to) the exact definition of 
25 the first peak, as well as the use of a mean sea surface height data set in place of 
the EGM08 geoid. Therefore, the comparisons done in this study are similar, but not 
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exact reproductions of methodologies. The purpose of the comparison is to highlight 
the physical basis between differences in the retracking methods. 

The mean difference in sea ice lead elevations retrieved by the waveform fitting 
method described in this study and the empirical tracker described in Giles et al. (2007) 
is 2.8cm (Giles et al., 2007 tracker-waveform fitting method), and the correlation is 
0.7. The most significant difference difference between the freeboard retrieval method 
of Laxon et al. (2013) and the waveform fitting method is the use of the 50% thresh- 
old tracker in the retrieval of sea ice floe elevations, which as illustrated in Sect. 3 
is expected to be biased high from theoretical arguments since the selected thresh- 
old should be closer to the waveform peak. Figure 11 shows the retrieved freeboard 
using the ELTF method, and Fig. 12 shows the difference with the waveform fitting 
method. The mean freeboard differences are 11.9 cm, 12.7 cm, and 11.5 cm for the 
March 2011-2013 periods. This corresponds to mean ice thickness differences of 
1.12m, 1.19m, and 1.08m using Eq. (17). The differences shown in Fig. 12 also show 
significant spatial and interannual differences between the methods. The mean free- 
board using the ELTF method for March 2011 is 31.3 cm, using Eq. (17) this corre- 
sponds to a minimum mean sea ice thickness of 2.9 m, which will be higher once one 
considers the contribution of snow. The waveform fitting method gives a minimum sea 
ice thickness of 1.8 m, which is much closer to the mean thickness of first year ice 
which is now the dominant ice type in the Arctic (Comiso, 2012). The h c snow speed 
of light correction was not applied in the comparison between threshold and waveform 
retrackers because it is equivalent for both data sets, but this will slightly increase the 
mean minimum thickness. Thus, the higher freeboard values retrieved by the threshold 
method are likely biased high, which is in agreement with the theoretical arguments 
presented in Sect. 3. 

5.2.2 Comparison to IceBridge data 

In order to compare the ELTF and waveform fitting methods, we now compare both 
methods to independent data collected from three measurement campaigns of the 
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Operation IceBridge mission. The mean IceBridge snow depth has been subtracted 
from the laser altimeter freeboard to determine the ice freeboard, and the data have 
then been gridded to the same 25 km polar stereographic grid as the CryoSat-2 data. 
A grid point is defined as containing valid data for comparison when there are greater 
5 than 200 IceBridge measurements and a valid gridded CryoSat-2 measurement. Since 
snow depth information is available from the IceBridge data set, we add the h c correc- 
tion factor to the CryoSat-2 retrieved freeboards and also estimate sea ice thickness 
using Eq. (17). 

Table 2 summarizes the comparison between the IceBridge observations, the wave- 
io form fitting method, and the ELTF method. For the waveform fitting method, the mean 
freeboard difference (CS2-lceBridge) ranges from 1-3 cm while the ice thickness dif- 
ference ranges from 11-23 cm. The slightly higher freeboard retrieved by CryoSat-2 
is consistent with the results of Armitage and Davidson, (2014) who estimate that the 
sea surface height will be biased low by ~ 2 cm due to off-nadir ranging to leads when 
is a minimum pulse peakiness of 0.18 is used as a threshold for lead classification. The 
mean freeboard difference (CS2-lceBridge) for the ELTF method ranges from 11.9- 
15.9 cm which corresponds to ice thickness differences of 11 2-1 49 cm, this is signifi- 
cantly higher than the waveform fitting method. As shown in Sect. 3, this is likely due 
to the choice of the 50% threshold which was shown to be too low in comparison to 
20 theoretical estimates which show the tracking point should be closer to the peak power 
value. Surface roughness and the finite sampling resolution of the radar also plays 
a role as well. We note that Laxon et al. (2013) did not add a correction for the speed 
of light within the snow pack and also subtracted a constant value from the sea surface 
elevation due to the use of different fitting models between open ocean and leads. In 
25 the waveform fitting retrieval scheme illustrated in this study, no such bias in the sea 
surface height needs to be removed because the same model is used to fit waveforms 
from open ocean, sea ice floes, and leads. The addition of the snow speed of light 
correction will also apply equally to each method. Thus, the waveform fitting method 
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gives a mean difference which compares much better to the IceBridge ice freeboard 
data using explicit geophysical arguments. 

The root-mean-square difference between the IceBridge data and CryoSat-2 free- 
board retrievals ranges from 7.4-11.1 cm for the waveform fitting method and is higher 
at 14.1-19.8 cm for the ELTF method. The mean estimated IceBridge freeboard uncer- 
tainty (taken from the data products using the method described in Kurtz et al., 2013) 
for the compared grid points is 5.9 cm, 7.6 cm, and 6.3 cm for the respective 2011-2013 
campaigns, the uncertainty in the sea surface height is due to a combination of instru- 
mental uncertainties and is also a function of distance to the nearest lead. For the 
waveform fitting method, the estimated uncertainty in the CryoSat-2 sea ice freeboard, 
(J CS 2 -fb> can be calculated as 

a cs2-fb = \J° diff - IceBridge C*8) 

where o diff is the observed standard deviation of differences and ff| Ce Bridge is the un ' 
certainty in the IceBridge ice freeboard. One complication with this estimate is that 
IceBridge ' s a set °f values, rather than a constant number. However, using the mean 
value of IceBridge f° r ea °h campaign gives an estimate for the CryoSat-2 freeboard 
uncertainty of 9.2 cm, 6.6 cm, and 3.8 cm for the respective 201 1-2013 campaigns. 

The correlation between the waveform fitting method and IceBridge data varies sub- 
stantially between the campaigns, but the correlation between the waveform fitting 
method and IceBridge data is higher for all campaigns than the ELTF method. While 
the low correlations to the IceBridge data may be some cause for concern, we note 
that the reasonable spatial distribution of sea ice freeboard shown in Fig. 9 and es- 
timated uncertainties for the IceBridge and CryoSat-2 freeboard retrievals place this 
into a context which can be understood. Given the uncertainties present in both the 
IceBridge and CryoSat-2 data and the small dynamic range of the freeboard values, it 
is possible a high correlation value can not be attained from the comparison. To test 
this hypothesis, a Monte Carlo simulation was conducted using the CryoSat-2 free- 
board values. Two sets of numbers were constructed by adding a random distribution 
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with zero mean and a standard deviation equivalent to the estimated mean uncertain- 
ties for the CryoSat-2 and IceBridge data. The correlation was computed for each set 
of numbers and the simulation was run 1 000 times. The simulation shows the expected 
correlations for two identical data sets with estimated uncertainties equivalent for the re- 
5 spective 201 1 , 2012, and 2013 campaigns are 0.46±0.05, 0.55±0.03, and 0.60±0.04. 
The March 2013 data set shows a correlation which is consistent with the estimated 
uncertainties for the data. However, it is not clear why only this campaign shows a lower 
RMS difference and a correlation which is in line with expectations. A possible explana- 
tion is that additional variability due to the use of non-temporally coincident data sets 
io added additional uncertainty to the estimates shown here, and that the lower bound 
uncertainty of 4cm is correct for the CryoSat-2 freeboard retrievals. A more detailed 
comparison between time coincident IceBridge data flights which underflew CryoSat-2 
will be the subject of a future study. 


6 Conclusions 

is A new method to fit CryoSat-2 level IB waveforms using an empirical model was de- 
veloped. This waveform fitting procedure was used to retrieve sea ice freeboard from 
CryoSat-2 over Arctic sea ice. Through comparison with Operation IceBridge data for 
the 2011-2013 campaigns, this study has shown that fitting of the CryoSat-2 level IB 
waveforms using a physical model can be used to obtain improved results over the em- 
20 pirical lead and threshold tracker (ELTF) methods which are similar to those used by 
Laxon et al. (2013). The ELTF method was found to have respective mean freeboard 
differences (CryoSat-2 - IceBridge) of 15.4 cm, 15.9 cm, and 11.9 cm and mean sea 
ice thickness differences of 144.2cm, 149.3cm, and 111.9cm. The mean freeboard 
differences for the waveform fitting method were 2.2 cm, 2.5 cm, and 1.1 cm, and the 
25 mean sea ice thickness differences were 20.6 cm, 23.3 cm, and 10.6 cm. The larger 
RMS and mean differences in the ELTF tracker method were found to be largely due to 
the choice of the 50% threshold, which was shown to be too low based on theoretical 
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modeling. The difference is also due to variations in surface roughness and the angular 
dependence of the backscattering coefficient. A bias of 1.9 cm was found in the wave- 
form fitting method freeboard retrievals compared to the IceBridge data, this bias is 
consistent with the estimated range bias due to off-nadir ranging of lead points shown 
5 by Armitage and Davidson (2014). A maximum correlation of 0.57 was found between 
the IceBridge freeboard and thickness data and the waveform fitting method, this cor- 
relation is consistent with an estimated uncertainty of 4 cm in the retrieved CryoSat-2 
freeboard for a 1 00 km gridded data point. 

Despite having a physical basis, and having a small bias compared to the airborne 
io observations, the model used to fit the CryoSat-2 waveforms is still essentially empiri- 
cal. In order to move towards a more physically exact model a number of points need 
to be taken into account, which are largely due to the considerable variability of sur- 
face types (and their associated backscattering properties) which can be found within 
the radar footprint: (1) it was assumed that the distribution of surface heights within 
is the footprint can be approximated with a Gaussian function, though it is possible that 
the presence of ridges will lead to a more skewed distribution. (2) The presence of 
ridges may also lead to an electromagnetic bias if the scattering from the ridge peaks 
is different than the surrounding ice, this is similar to a known phenomenon which has 
been observed in open ocean returns wherein wave troughs have a higher backscatter 
20 than wave crests (Yaplee et al., 1971). (3) The model assumes the antenna boresight 
is always at nadir and the surface normal is parallel to the nadir direction. However, re- 
cently discovered pitch and roll biases within the CryoSat-2 data mean that the antenna 
boresight is slightly off-nadir which should be taken into account. (4) Mixed returns 
containing more than one surface type are not dealt with in the model, and the fitting 
25 procedure only works when the mixture of surface types allows for sufficient separation 
between peaks to fit the return. The usage of additional statistics such as goodness 
of fit estimates may be used to further reduce errors caused by mixed returns in fu- 
ture studies. (5) Surface scattering from the snow-air interface, and volume scattering 
within the snow pack and sea ice were not considered. In particular, this may cause 
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additional uncertainty if the snow cover contains dense layers or if the roughness of 
the snow-air interface is much lower than that of the snow-ice interface. 

Further improvements to the retrieval of surface elevation, sea ice freeboard, and 
sea ice thickness can also be done which should reduce the uncertainty in the mea- 
surements. Modeling of the phase and further analysis of the SARIn data areas in the 
Arctic Ocean may lead to ways to identify off-nadir sea ice leads and reduce the ob- 
served ~2cm bias which this causes in the sea surface height and freeboard data 
sets. In a similar manner to what was described in Laxon et al. (2013), the retrieval of 
a large volume of sea surface height estimates will allow for the construction of a high 
resolution mean sea surface height data set which can be used to reduce geoid errors 
which are known to be prevalent in the Arctic Ocean (McAdoo et al., 2013). Lastly, an 
evaluation of the IceBridge snow depth measurements needs to be done to improve 
basin-wide snow depth on sea ice estimates. This has been done for a single season 
of data (Kurtz and Farrell, 2011) compared to the snow depth climatology of Warren 
et al. (1 999), and for passive microwave retrievals of snow depth on sea ice for first year 
ice (Brucker and Markus, 2013). The focus of a future study will be to utilize existing 
observations to improve estimates of snow depth on sea ice to be used in the retrieval 
of sea ice thickness from the CryoSat-2 time series. 

Overall, this study has further demonstrated the capabilities of CryoSat-2 for the 
retrieval of sea ice freeboard and thickness. The advantage of the retrieval processes 
used in this study is that they are compatible with the laser altimetry record and show 
that the two records can be reconciled to produce a more complete time series of sea 
ice volume change. This has distinct advantages for the expected launch of the ICESat- 
2 laser altimeter mission in 2017. The lifetime of CryoSat-2 is expected to overlap with 
the ICESat-2 mission, as is the new Sentinel-3 radar altimeter mission. The combined 
satellite radar and laser altimetry data provided by these missions will thus provide 
unmatched information on the state of the Arctic sea ice cover. 
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Table 1. Summary of parameters and symbols used in the CryoSat-2 model. 


Symbol in text 

Parameter 

Average value 

A 

Center wavelength 

0.0221 m 

G 0 

One-way antenna gain 

42 dB 

D 0 

One-way gain of synthetic beam 

30.6 dB 

c 

Speed of light in vacuum 

299 792 458 ms -1 

o°(0°) 

Nadir backscatter coefficient 


h 

Satellite altitude 

725 km 

n 

Geometric factor 

1.113 

N b 

Number of synthetic beams 

64 

T 

Echo delay time 


tk 

Look angle of synthetic beam k from nadir 


H 

Heaviside step function 


Yi 

Elliptical antenna pattern term 1 

6767.6 

v 2 

Elliptical antenna pattern term 2 

664.06 

a 

Angular backscattering efficiency 


^0 

Carrier wavenumber 

284.307 rrf 1 


Satellite velocity 

7435 ms -1 

a 

Standard deviation of surface height 


fi w 

Received bandwidth 

320 MHz 


g 

C/5* 

o 
c 
c n 
co 
O' 


"0 


Q) 

"O 

CD 


CryoSat-2 sea ice 
freeboard 


a 


< / > 
O 

c 
c n 

(D 


TJ 

o 

■o 

CD 


g 

C/)’ 

o 

c 

c/) 

C/) 

o' 


TJ 

0) 

"O 

CD 


D 

c/)’ 

o 

c: 

(j) 

(j) 

o’ 

Z3 

"D 

0) 

"O 

CD 


N. T. Kurtz et al. 


Title Page 


Abstract H Introduction 


Conclusions H References 


Figures 



Full Screen / Esc 


Printer-friendly Version 
Interactive Discussion 



755 

















TCD 

8, 721-768, 2014 


g 

C/5* 

o 
c 
c n 
co 
O' 


"0 


03 

"O 

CD 


CryoSat-2 sea ice 
freeboard 


g 

(/)’ 

O 

c 
c n 
V) 

O' 


Table 2. CryoSat-2 freeboard and thickness retrievals compared to IceBridge airborne data. _g 
Values pertaining to sea ice thickness are in parentheses. 



Mar 2011 

Mar 2012 

Mar 2013 

Number of grid points 

223 


391 


217 



Waveform fitting 

ELTF 

Waveform fitting 

ELTF 

Waveform fitting 

ELTF 

Mean difference (CS2 - IceBridge) 

2.2 cm (20.6 cm) 

15.4cm (144.2cm) 

2.5 cm (23.3 cm) 

15.9cm (149.3cm) 

1.2 cm (11.4 cm) 

11.9 cm (111.9 cm) 

RMS difference 

11.1 cm 

19.7 cm 

10.4 cm 

19.8 cm 

7.4 cm 

14.1 cm 

Correlation 

0.02 

-0.12 

0.25 

0.11 

0.57 

0.55 
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Fig. 1. Example CryoSat-2 waveforms, (a) Example waveform for a sea ice floe. The waveform 
contains an off-nadir reflection from a surface with high backscatter which results in a secondary 
peak, (b) Example specular waveform from a lead. 
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Fig. 2. Maps of IceBridge ice freeboard data for 16-28 March 2011, 14 March-10 April 2012, 
and 20-27 March 2013. 
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Echo delav fnsl 



Fig. 3. Simulated CryoSat-2 echoes over sea ice leads for a = 0.02 m and (a) the typical range 
of a observed over sea ice floes, (b) a zoomed in plot showing the behavior near the echo 
delay time of 0. 
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Fig. 4. Simulated CryoSat-2 echoes over sea ice floes for (a) the typical range of delay times 
provided in the level IB product over sea ice, (b) a zoomed in plot showing the behavior near 
the echo delay time of 0. Solid lines correspond to a = 10 3 and dashed lines correspond to 
a= 10 5 . 
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Fig. 5. Ratio of the average power of points located within 10 ns following the point of peak 
power, to the peak power. These results are taken directly from the physical model in Sect. 3. 
This is used to provide an initial guess in the fitting of waveforms over leads. The x axis is 
a logarithmic scale to better show the variation over the large dynamic range of a. 
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Fig. 6. Example CryoSat-2 waveform fits for sea ice leads with increasing values for a. The g 
fitted waveform at the instrument sampling resolution is shown in blue, and the CryoSat-2 data 
points are represented by black dots, (a) a = 5.1 xIO 5 , o - 0.01 m. (b) a = 2.7x 10 6 , o = 0.01 m. 

(c) a - 2.8 X 1 0 7 , a = 0.007 m. (d) a = 5.0 x 1 0 7 , <r = 0.005 m. 
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Fig. 7. Ratio of the average power of points located within 90 to 120 ns following the point of 
peak power, to the peak power. These results are taken directly from the physical model in 
Sect. 3. This is used to provide an initial guess in the fitting of waveforms of sea ice floes. The 
x axis is a logarithmic scale to better show the variation over the large dynamic range of a. 
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Fig. 8. Example CryoSat-2 waveform fits for sea ice floes. The fitted waveform is shown in blue 
and the CryoSat-2 data are represented by black dots, (a) and (b) are fits of waveforms which 
demonstrate the good agreement between the observations and the model, (c) and (d) are 
fits of waveforms containing multiple peaks in the trailing edge due to the presence of strong 
off-nadir reflections from smooth ice and/or leads. 
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Fig. 9. Gridded freeboard retrievals from the CryoSat-2 waveform fitting method and their dis- 
tributions. 
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Fig. 10. Map of the mean gridded surface roughness, a, from CryoSat-2 excluding sea ice 
lead points. Surface roughness from low ice concentration areas near the sea ice edge have 
also been included, these areas have a high surface roughness due to the presence of ocean 
waves. 
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Fig. 11. Maps of the monthly mean gridded sea ice freeboard, fb radar , and their probability 
distributions from CryoSat-2 using a 50% threshold tracker for sea ice floes similar to Laxon 2 
et al. (2013) and an empirical retracker for sea ice leads used by Giles et al. (2007). 
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Fig. 12. Freeboard difference between the waveform fitting method and the ELTF method. 
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